Brief description of the script
This R markdown document reads, summarizes and plots data for: Paixao et al. 2021. Functional analysis of the Middle Paleolithic Ground Stone Tools from Unit V of Nesher Ramla (Central Levante): the application of a multi-scale use-wear approach. Quaternary International
The document contains includes plots of the quantitative surface texture analysis, using Confocal microcopy.
This R project and respective scripts follow the procedures described by Marwick et al. 2017.
To compile this markdown document do not delete or move files from their original folders. Please note that most of the tables and figures in this file do not match the numbering in the PhD dissertation manuscript.
For any questions, comments and inputs, please contact:
Eduardo Paixão, paixao@rgzm.de
Imported files are in: ‘../analysis/raw_data’
Figures are saved in: ‘../analysis/plots’
Tables are saved in: ‘../analysis/derived_data’
# Load required libraries
library(tidyverse)
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.0 ──
## ✓ ggplot2 3.3.3 ✓ purrr 0.3.4
## ✓ tibble 3.1.0 ✓ dplyr 1.0.5
## ✓ tidyr 1.1.3 ✓ stringr 1.4.0
## ✓ readr 1.4.0 ✓ forcats 0.5.1
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
library(utils)
library(knitr)
library(janitor)
##
## Attaching package: 'janitor'
## The following objects are masked from 'package:stats':
##
## chisq.test, fisher.test
library(kableExtra)
##
## Attaching package: 'kableExtra'
## The following object is masked from 'package:dplyr':
##
## group_rows
library(GGally)
## Registered S3 method overwritten by 'GGally':
## method from
## +.gg ggplot2
library(doBy)
##
## Attaching package: 'doBy'
## The following object is masked from 'package:dplyr':
##
## order_by
library(ggpubr)
library(tools)
# See your WD and update the following paths
# getwd()
# Load data from .csv
confocaldataarch <- read.delim("../raw_data/confocalarch/confocaldataarch.csv", header = T, ";")
data_file <- list.files("../raw_data/confocalarch", pattern = "\\.csv$", full.names = TRUE)
md5_in <- md5sum(data_file)
info_in <- data.frame(file = basename(names(md5_in)), checksum = md5_in, row.names = NULL)
# compute descriptive statistics
nminmaxmeanmedsd <- function(x){
y <- x[!is.na(x)]
n_test <- length(y)
min_test <- min(y)
max_test <- max(y)
mean_test <- mean(y)
med_test <- median(y)
sd_test <- sd(y)
out <- c(n_test, min_test, max_test, mean_test, med_test, sd_test)
names(out) <- c("n", "min", "max", "mean", "median", "sd")
return(out)
}
num.var <- 21:length(confocaldataarch)
confostatsarch <- summaryBy(.~sample + workedmaterial, data=confocaldataarch[c("sample", "workedmaterial", names(confocaldataarch)[num.var])], FUN=nminmaxmeanmedsd)
write_csv(confostatsarch, "../derived_data/confocalstats_arch.csv")
# Sample main dataset
# Only archaeological tools
confoarch <- filter(confocaldataarch, sample == "archaeological")
# Loop for plotting all surface texture parameters
for (i in num.var) cat("[",i,"] ", names(confoarch)[i], "\n", sep = "")
## [21] Sq
## [22] Ssk
## [23] Sku
## [24] Sp
## [25] Sv
## [26] Sz
## [27] Sa
## [28] Smr
## [29] Smc
## [30] Sxp
## [31] Sal
## [32] Str
## [33] Std
## [34] Sdq
## [35] Sdr
## [36] VM
## [37] Vv
## [38] Vmp
## [39] Vmc
## [40] Vvc
## [41] Vvv..p...80.00..
## [42] Vvv
## [43] Mean.depth.of.furrows
## [44] Mean.density.of.furrows
## [45] First.direction
## [46] Second.direction
## [47] Third.direction
## [48] Isotropy
## [49] Lengthscale.anisotropy.Sfrax.epLsar
## [50] Length.scale.anisotropy..NewEplsar.
## [51] Fractal.complexity.Asfc
## [52] Smfc
## [53] HAsfc9
## [54] HAsfc81
for (i in num.var) {
p <- ggplot(data = confocaldataarch, aes_string(x = "workedmaterial", y = names(confoarch)[i])) +
geom_boxplot() +
# geom_line(aes(group = motion)) +
theme_classic() +
# facet_wrap(~ sample) +
labs(x = "material", y = gsub("\\.", " ", names(confoarch)[i])) +
scale_colour_hue(h = c(25,225))
print(p)
# saves the plots
file_out <- paste0(file_path_sans_ext(info_in[["file"]]), "_plot_",
names(confoarch)[i], ".pdf")
ggsave(filename = file_out, plot = p, path = "../plots", device = "pdf", width = 26,
height = 21, units = "cm" )
}
# Comparing archaeological and experimental
# Loop for plotting all surface texture parameters
for (i in num.var) cat("[",i,"] ", names(confocaldataarch)[i], "\n", sep = "")
## [21] Sq
## [22] Ssk
## [23] Sku
## [24] Sp
## [25] Sv
## [26] Sz
## [27] Sa
## [28] Smr
## [29] Smc
## [30] Sxp
## [31] Sal
## [32] Str
## [33] Std
## [34] Sdq
## [35] Sdr
## [36] VM
## [37] Vv
## [38] Vmp
## [39] Vmc
## [40] Vvc
## [41] Vvv..p...80.00..
## [42] Vvv
## [43] Mean.depth.of.furrows
## [44] Mean.density.of.furrows
## [45] First.direction
## [46] Second.direction
## [47] Third.direction
## [48] Isotropy
## [49] Lengthscale.anisotropy.Sfrax.epLsar
## [50] Length.scale.anisotropy..NewEplsar.
## [51] Fractal.complexity.Asfc
## [52] Smfc
## [53] HAsfc9
## [54] HAsfc81
for (i in num.var) {
p <- ggplot(data = confocaldataarch, aes_string(x = "workedmaterial", y = names(confocaldataarch)[i],
colour = "sample")) +
geom_boxplot() +
# geom_line(aes(group = motion)) +
theme_classic() +
labs(colour = "sample") +
# facet_wrap(~ sample) +
labs(x = "material", y = gsub("\\.", " ", names(confocaldataarch)[i])) +
scale_colour_hue(h = c(25,225), limits = levels(confocaldataarch[["sample"]]))
print(p)
# saves the plots
file_out <- paste0(file_path_sans_ext(info_in[["file"]]), "_plot_",
names(confocaldataarch)[i], ".pdf")
ggsave(filename = file_out, plot = p, path = "../plots", device = "pdf", width = 26,
height = 21, units = "cm" )
}
# Only archaeological
# Sa vs. Sq
Sa_Sq <- ggplot(data = confoarch) +
geom_point(mapping = aes(x = Sa, y = Sq, colour = workedmaterial)) +
theme_classic() +
labs(colour = "workedmaterial") +
scale_colour_hue(h = c(25, 230))
print(Sa_Sq)
file_out <- paste0(file_path_sans_ext(info_in[["file"]]), "_scatterplot_Sa-Sq", ".pdf")
ggsave(filename = file_out, plot = Sa_Sq, path = "../plots", device = "pdf")
## Saving 7 x 5 in image
# epLsar vs. Asfc
ep_As <- ggplot(data = confoarch) +
geom_point(mapping = aes(x = Fractal.complexity.Asfc, y = Lengthscale.anisotropy.Sfrax.epLsar, colour = workedmaterial)) +
theme_classic() +
labs(colour = "workedmaterial") +
scale_colour_hue(h = c(25, 230))
print(ep_As)
file_out <- paste0(file_path_sans_ext(info_in[["file"]]), "_scatterplot_Asfc-epLsar", ".pdf")
ggsave(filename = file_out, plot = ep_As, path = "../plots", device = "pdf")
## Saving 7 x 5 in image
# Sq vs. Vmc
Sq_Vmc <- ggplot(data = confoarch) +
geom_point(mapping = aes(x = Sq, y = Vmc, colour = workedmaterial)) +
theme_classic() +
labs(colour = "workedmaterial") +
scale_colour_hue(h = c(25, 230))
print(Sq_Vmc)
file_out <- paste0(file_path_sans_ext(info_in[["file"]]), "_scatterplot_Sq-Vmc", ".pdf")
ggsave(filename = file_out, plot = Sq_Vmc, path = "../plots", device = "pdf")
## Saving 7 x 5 in image
# Mean depth of furrows vs. mean density of furrows
furrows <- ggplot(data = confoarch) +
geom_point(mapping = aes(x = Mean.depth.of.furrows, y = Mean.density.of.furrows,
colour = workedmaterial)) +
theme_classic() +
labs(colour = "workedmaterial", x = "Mean depth of furrows", y = "Mean density of furrows") +
scale_colour_hue(h = c(25, 230))
print(furrows)
file_out <- paste0(file_path_sans_ext(info_in[["file"]]), "_scatterplot_furrows", ".pdf")
ggsave(filename = file_out, plot = furrows, path = "../plots", device = "pdf")
## Saving 7 x 5 in image
# combine all in a single image
ggarrange(Sa_Sq, Sq_Vmc, furrows, ep_As, common.legend = TRUE, legend = "bottom")
ggsave("../plots/scatterplots.png")
## Saving 7 x 5 in image
# Comparing archaeological and experimental
# Sa vs. Sq
Sa_Sq <- ggplot(data = confocaldataarch) +
geom_point(mapping = aes(x = Sa, y = Sq, colour = workedmaterial)) +
theme_classic() +
labs(colour = "workedmaterial") +
facet_wrap(~ sample) +
scale_colour_hue(h = c(25, 230))
print(Sa_Sq)
file_out <- paste0(file_path_sans_ext(info_in[["file"]]), "_scatterplot_Sa-Sq", ".pdf")
ggsave(filename = file_out, plot = Sa_Sq, path = "../plots", device = "pdf")
## Saving 7 x 5 in image
# epLsar vs. Asfc
ep_As <- ggplot(data = confocaldataarch) +
geom_point(mapping = aes(x = Fractal.complexity.Asfc, y = Lengthscale.anisotropy.Sfrax.epLsar, colour = workedmaterial)) +
theme_classic() +
labs(colour = "workedmaterial") +
facet_wrap(~ sample) +
scale_colour_hue(h = c(25, 230))
print(ep_As)
file_out <- paste0(file_path_sans_ext(info_in[["file"]]), "_scatterplot_Asfc-epLsar", ".pdf")
ggsave(filename = file_out, plot = ep_As, path = "../plots", device = "pdf")
## Saving 7 x 5 in image
# Sq vs. Vmc
Sq_Vmc <- ggplot(data = confocaldataarch) +
geom_point(mapping = aes(x = Sq, y = Vmc, colour = workedmaterial)) +
theme_classic() +
labs(colour = "workedmaterial") +
facet_wrap(~ sample) +
scale_colour_hue(h = c(25, 230))
print(Sq_Vmc)
file_out <- paste0(file_path_sans_ext(info_in[["file"]]), "_scatterplot_Sq-Vmc", ".pdf")
ggsave(filename = file_out, plot = Sq_Vmc, path = "../plots", device = "pdf")
## Saving 7 x 5 in image
# Mean depth of furrows vs. mean density of furrows
furrows <- ggplot(data = confocaldataarch) +
geom_point(mapping = aes(x = Mean.depth.of.furrows, y = Mean.density.of.furrows,
colour = workedmaterial)) +
theme_classic() +
labs(colour = "workedmaterial", x = "Mean depth of furrows", y = "Mean density of furrows") +
facet_wrap(~ sample) +
scale_colour_hue(h = c(25, 230))
print(furrows)
file_out <- paste0(file_path_sans_ext(info_in[["file"]]), "_scatterplot_furrows", ".pdf")
ggsave(filename = file_out, plot = furrows, path = "../plots", device = "pdf")
## Saving 7 x 5 in image
# combine all in a single image
ggarrange(Sa_Sq, Sq_Vmc, furrows, ep_As, common.legend = TRUE, legend = "bottom")
ggsave("../plots/scatterplots.png")
## Saving 7 x 5 in image
data(confocaldataarch, package = "reshape")
## Warning in data(confocaldataarch, package = "reshape"): data set
## 'confocaldataarch' not found
data(confoarch, package = "reshape")
## Warning in data(confoarch, package = "reshape"): data set 'confoarch' not found
# Only archaeological
# Height parameters
ggpairs(data=confoarch,
columns = c(21:26),
cardinality_threshold = 30,
mapping = ggplot2::aes(color = workedmaterial),
lower = list(continuous = wrap("points", alpha = 0.5, size = 1)),
upper = list(continuous = "blank"),
legend = c(2,1)
) +
theme(legend.position = "right") +
labs(fill = "Micro polish type")
ggsave("../plots/confocalarcharea_matrix.png")
## Saving 7 x 5 in image
# Volume parameters
ggpairs(data=confoarch,
columns = c(36:41),
cardinality_threshold = 30,
mapping = ggplot2::aes(color = workedmaterial),
lower = list(continuous = wrap("points", alpha = 0.5, size = 1)),
upper = list(continuous = "blank"),
legend = c(2,1)
) +
theme(legend.position = "right") +
labs(fill = "Micro polish type")
ggsave("../plots/confocalarchvolume_matrix.png")
## Saving 7 x 5 in image
# Comparing archaeological and experimental
# Height parameters
ggpairs(data=confocaldataarch,
columns = c(21:26),
cardinality_threshold = 30,
mapping = ggplot2::aes(color = workedmaterial, shape = sample),
lower = list(continuous = wrap("points", alpha = 0.5, size = 1)),
upper = list(continuous = "blank"),
legend = c(2,1)
) +
theme(legend.position = "right") +
labs(fill = "Micro polish type")
ggsave("../plots/confocalarcharea_matrix.png")
## Saving 7 x 5 in image
# Volume parameters
ggpairs(data=confocaldataarch,
columns = c(36:41),
cardinality_threshold = 30,
mapping = ggplot2::aes(color = workedmaterial, shape = sample),
lower = list(continuous = wrap("points", alpha = 0.5, size = 1)),
upper = list(continuous = "blank"),
legend = c(2,1)
) +
theme(legend.position = "right") +
labs(fill = "Micro polish type")
ggsave("../plots/confocalarchvolume_matrix.png")
## Saving 7 x 5 in image
# select parameter from dataset
# Only archaeological
confostatsarch2 <- filter(confostatsarch, sample == "archaeological")
heightconfostats <- select(confostatsarch2,workedmaterial, Sq.mean,Ssk.mean,Sku.mean,Sp.mean,Sv.mean,Sz.mean,Sa.mean)
p1 <- ggplot(heightconfostats, aes(x=workedmaterial, y=Sq.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p2 <- ggplot(heightconfostats, aes(x=workedmaterial, y=Ssk.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p3 <- ggplot(heightconfostats, aes(x=workedmaterial, y=Sku.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p4 <- ggplot(heightconfostats, aes(x=workedmaterial, y=Sp.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p5 <- ggplot(heightconfostats, aes(x=workedmaterial, y=Sv.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p6 <- ggplot(heightconfostats, aes(x=workedmaterial, y=Sz.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p7 <- ggplot(heightconfostats, aes(x=workedmaterial, y=Sa.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
ggarrange(p1, p2, p3, p4, p5, p6, p7, common.legend = TRUE, font.label = list(size=8), legend="bottom")
ggsave("../plots/confostatsarcharea_boxplots.png")
## Saving 7 x 5 in image
# Now Volume parameters
volumeconfostats <- select(confostatsarch,sample,workedmaterial, VM.mean,Vv.mean,Vmp.mean,Vmc.mean,Vvc.mean,Vvv.mean)
p8 <- ggplot(volumeconfostats, aes(x=workedmaterial, y=VM.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p9 <- ggplot(volumeconfostats, aes(x=workedmaterial, y=Vv.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p10 <- ggplot(volumeconfostats, aes(x=workedmaterial, y=Vmp.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p11 <- ggplot(volumeconfostats, aes(x=workedmaterial, y=Vmc.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p12 <- ggplot(volumeconfostats, aes(x=workedmaterial, y=Vvc.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p13 <- ggplot(volumeconfostats, aes(x=workedmaterial, y=Vvv.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
ggarrange(p8, p9, p10, p11, p12, p13, common.legend = TRUE, font.label = list(size=8), legend="bottom")
ggsave("../plots/confostatarchvolume_boxplots.png")
## Saving 7 x 5 in image
# Comparing archaeological and experimental
# first Height parameters
heightconfostats <- select(confostatsarch,sample,workedmaterial, Sq.mean,Ssk.mean,Sku.mean,Sp.mean,Sv.mean,Sz.mean,Sa.mean)
p1 <- ggplot(heightconfostats, aes(x=workedmaterial, y=Sq.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p2 <- ggplot(heightconfostats, aes(x=workedmaterial, y=Ssk.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p3 <- ggplot(heightconfostats, aes(x=workedmaterial, y=Sku.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p4 <- ggplot(heightconfostats, aes(x=workedmaterial, y=Sp.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p5 <- ggplot(heightconfostats, aes(x=workedmaterial, y=Sv.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p6 <- ggplot(heightconfostats, aes(x=workedmaterial, y=Sz.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p7 <- ggplot(heightconfostats, aes(x=workedmaterial, y=Sa.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
ggarrange(p1, p2, p3, p4, p5, p6, p7, common.legend = TRUE, font.label = list(size=8), legend="bottom")
ggsave("../plots/confostatsarcharea_boxplots.png")
## Saving 7 x 5 in image
# Now Volume parameters
volumeconfostats <- select(confostatsarch,sample,workedmaterial, VM.mean,Vv.mean,Vmp.mean,Vmc.mean,Vvc.mean,Vvv.mean)
p8 <- ggplot(volumeconfostats, aes(x=workedmaterial, y=VM.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p9 <- ggplot(volumeconfostats, aes(x=workedmaterial, y=Vv.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p10 <- ggplot(volumeconfostats, aes(x=workedmaterial, y=Vmp.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p11 <- ggplot(volumeconfostats, aes(x=workedmaterial, y=Vmc.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p12 <- ggplot(volumeconfostats, aes(x=workedmaterial, y=Vvc.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
p13 <- ggplot(volumeconfostats, aes(x=workedmaterial, y=Vvv.mean, colour=workedmaterial)) +
geom_boxplot() +
labs(x="", colour="Micro polish")
ggarrange(p8, p9, p10, p11, p12, p13, common.legend = TRUE, font.label = list(size=8), legend="bottom")
ggsave("../plots/confostatarchvolume_boxplots.png")
## Saving 7 x 5 in image
sessionInfo()
## R version 4.0.4 (2021-02-15)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Catalina 10.15.7
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRblas.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] tools stats graphics grDevices utils datasets methods
## [8] base
##
## other attached packages:
## [1] ggpubr_0.4.0 doBy_4.6.9 GGally_2.1.1 kableExtra_1.3.4
## [5] janitor_2.1.0 knitr_1.31 forcats_0.5.1 stringr_1.4.0
## [9] dplyr_1.0.5 purrr_0.3.4 readr_1.4.0 tidyr_1.1.3
## [13] tibble_3.1.0 ggplot2_3.3.3 tidyverse_1.3.0
##
## loaded via a namespace (and not attached):
## [1] fs_1.5.0 lubridate_1.7.10 webshot_0.5.2 RColorBrewer_1.1-2
## [5] httr_1.4.2 Deriv_4.1.3 backports_1.2.1 bslib_0.2.4
## [9] utf8_1.2.1 R6_2.5.0 DBI_1.1.1 colorspace_2.0-0
## [13] withr_2.4.1 tidyselect_1.1.0 gridExtra_2.3 curl_4.3
## [17] compiler_4.0.4 cli_2.3.1 rvest_1.0.0 xml2_1.3.2
## [21] labeling_0.4.2 sass_0.3.1 scales_1.1.1 systemfonts_1.0.1
## [25] digest_0.6.27 foreign_0.8-81 rmarkdown_2.7 svglite_2.0.0
## [29] rio_0.5.26 pkgconfig_2.0.3 htmltools_0.5.1.1 dbplyr_2.1.0
## [33] highr_0.8 rlang_0.4.10 readxl_1.3.1 rstudioapi_0.13
## [37] jquerylib_0.1.3 generics_0.1.0 farver_2.1.0 jsonlite_1.7.2
## [41] zip_2.1.1 car_3.0-10 magrittr_2.0.1 Matrix_1.3-2
## [45] Rcpp_1.0.6 munsell_0.5.0 fansi_0.4.2 abind_1.4-5
## [49] lifecycle_1.0.0 stringi_1.5.3 yaml_2.2.1 snakecase_0.11.0
## [53] carData_3.0-4 MASS_7.3-53.1 plyr_1.8.6 grid_4.0.4
## [57] crayon_1.4.1 lattice_0.20-41 haven_2.3.1 cowplot_1.1.1
## [61] hms_1.0.0 pillar_1.5.1 ggsignif_0.6.1 reprex_1.0.0
## [65] glue_1.4.2 evaluate_0.14 data.table_1.14.0 modelr_0.1.8
## [69] vctrs_0.3.6 cellranger_1.1.0 gtable_0.3.0 reshape_0.8.8
## [73] assertthat_0.2.1 xfun_0.22 openxlsx_4.2.3 broom_0.7.5
## [77] rstatix_0.7.0 viridisLite_0.3.0 ellipsis_0.3.1